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Abstract 

We introduce a family of models for heat conduction with and without momentum 
conservation. They are analytically solvable in the high temperature limit and can 
also be efficiently simulated. In all cases Fourier Law is verified in one dimension. 
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1 Introduction 



When there is a temperature difference between the boundaries of a material, heat is 
transported from the hottest to the coldest side. The phenomenological law governing 
this process has been known for a long time: the Fourier law J = kVT states the 
proportionality of the heat flux J (the amount of heat transported through the unit 
surface in unit time) to thermal gradient VT (the spatial derivative of the temperature 
field). The proportionality constant k is called thermal conductivity coefficient. Almost 
two centuries after Fourier's law was discovered, its microscopic derivation is still an 
open problem from a fundamental point of view. At stake is not only a question of 
mathematical rigor: spatial constraints in certain experimental systems can significantly 
alter the transport properties in ways that are not yet fully understood. Systems with 
dimension d < 2 can have a thermal conductivity that can even become anomalous (i.e. 
diverging in the thermodynamical limit), implying a breakdown of the usual Fourier law. 
In the last few years such systems became experimentally realizable, and this has spurred a 
renewed interest in the theoretical basis among the physics community (for recent reviews 
see |UE!). 

An important family of models which has been extensively studied consist of regular 
periodic d dimensional lattices of N = L d point-like atoms interacting with their neighbors 
through non-linear forces, as in the Fermi-Pasta-Ulam (FPU) chain. For the linear chain 
(d = 1) of oscillators it has been rigorously established in the sixties |3] that the heat 
conductivity diverges in the thermodynamic limit like k ~ L. This result was later 
generalized to the higher dimensional case in 

Indeed, the dynamics of such a system can be decomposed into the evolution of non- 
interacting waves (normal modes) which do not exchange energy amongst each other, 
thus giving rise to ballistic transport. It is nowadays clear that in order to have a fi- 
nite thermal coefficient the model must posess an efficient scattering mechanism between 
phonons. However, numerical studies indicate that also in the presence of anharmonicity 
heat conductivity is still anomalous, with a power law divergence k ~ L a , where a ~ 0.37 
jj] for a linear chain. In recent years, an increasingly large amount of (very often con- 
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flicting) numerical results have appeared in the literature. Although some progress have 
been achieved, many puzzles remain. At present, the picture we have for Id is as follows: 

• Anomalous conductivity occurs in cases in which momentum is conserved, i.e. at 
least one acoustic phonon branch is present in the harmonic limit |5|IF)|IT|IH|ll(}|lllj. 
The reason of the anomalous behavior has been traced back to the low-frequency 
Fourier modes, which are only weakly damped by the interaction with other modes 
so they behave as very efficient energy carriers for the system. The condition of 
momentum-conservation is however not sufficient: the 'coupled rotor model' seems 
to have normal conductivity O EH, the reason being attributed to phase jumps 
between barriers of the periodic potential which act as scatters even for the long 
wave-length modes. 

Here a short discussion is in order: if a system is composed of free particles, mo- 
mentum conservation is associated with the symmetry with respect to simultaneous 
translation of all particles. In a system such as the FPU or the coupled rotator 
system, we have particles with coordinates Xi - a deviation from the lattice position 
in the former, and an angle in the latter case. By translation invariance we may 
mean the symmetry 'in the coordinates' i — *■ i + 1, or the symmetry 'in the fields' 
Xi — > Xi + 5. The issue is somewhat mixed up by the fact that in the FPU system we 
tend to picture the Xj along the line joining the baths, while for the coupled rotors 
we think of the X{ as being transversal. In fact, as already pointed out in [33], the 
difference between 'longitudinal' and 'transversal' only arises if the thermal baths 
interact with all the particles which reach a certain position (so that there is contact 
with the bath or not depending on the actual value of Xi), while it is irrelevant if 
the bath is connected to some fixed particles at the end, whatever their positions 
Xi might be. One could conjecture that that the transverse nature of the Xj in the 
rotor case was responsible for the finite conductivity, although this would suggest 
that the FPU chain with baths coupled to fixed particles at the ends should have 
normal conductivity too, contrary to observation. 
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• Finite thermal conductivity has been obtained for special models with interaction 
with a substrate which violates momentum conservation (the so-called Ding-a-ling 
[Tf] and Ding-dong 15j models or various Klein-Gordon chains I17j). 

• A Hamiltonian system for which a macroscopic transport law has been rigorously 
derived is a gas of noninteracting particles moving among a fixed array of periodic 
convex scatters (periodic Lorentz gas or Sinai billiard) [T% | ITT H I2T)]. More recently, 
a modified Lorentz gas with fixed freely-rotating circular scatters has also been 
considered |2lj. In a further development, Eckmann and Young introduced a 
related model (having also normal heat conduction) that can be solved analytically 
in a limit. 

• The consequences of the underlying microscopic dynamics are very controversial: 
while it was believed for a long time that dynamical chaos and full hyperbolicity 
is a necessary and sufficient condition in order to have finite conductivity, recent 
examples show that: i) deterministic diffusion and normal heat transport can take 
place in systems with linear instability [2211231121]; h) there are systems with positive 
Lyapunov exponents (FPU model) where the heat conduction does not obey Fourier 
law. 

• The role of disorder (in the form of random masses of the atoms) has also been 
considered. For the linear oscillators it has been rigorously shown that disorder is 
not enough to reproduce correctly Fourier law (in fact one finds k ~ L 1 / 2 [2I3l2nil2Zl)- 
For non-linear oscillators the anomalous behavior has been confirmed numerically 

RSI- 

• Stochastic models have also been studied fZ§\ E?Uj and a finite thermal conductivity 
has been found. The large deviation properties of this kind of models have also been 
investigated recently [32] • While these models can be useful to obtain exact results 
under the assumption of suitable Markovian dynamics, one would like to obtain an 
explanation of Fourier law starting from a pure deterministic description. 
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In the presence of anomalous heat conductivity there exist mainly two general theoretical 
predictions: 

• A self consistent 'mode-coupling' theory which predicts an exponent a = 2/5 [33| 12] 

• Renormalization group analysis of a set of hydrodynamic equations for a Id fluid 
which yields a = 1/3 [2]. It has been suggested that a crossover might exist 
between these two behaviours, see [35] . 

Both of them apply to momentum conserving systems. It is worth mentioning that present 
numerical simulations are not fully in agreement with either mode-coupling or liquid the- 
ory predictions. Some authors suggested that the discrepancy might arise because inter- 
mediate sizes and timescales currently available in simulations are not able to probe the 
real behavior in the thermodynamic limit. In this respect, the problem must be admit- 
tedly seen as an open question. The situation is even less clear in 2c? lattices. Theoretical 
arguments support a weak logarithmic divergence for the thermal conductivity. On the 
other hand numerical simulations gave very conflicting results: logarithmic divergence 
[35] . power law anomalous behavior [37] and even finite thermal conductivity |38j . 

To understand some of the basic features of the energy transport in l-d systems, we 
develop in this paper a simple model, which can be cast also in the form of a set of discrete 
coupled maps such that in a limit the continuum version is reobtained. It can be solved 
analytically for high energies, and it is suitable at all energies for numerical simulations 
because it is much less computer-time demanding compared to a continuous flux. Because 
the approach to the high-energy value of the conductivity can be seen in the simulations, 
the discussion is not plagued with doubts as to whether the results are preasymptotic. 
We have found no differences between the case in which momentum is conserved and the 
case in which it is not. 
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2 Hamiltonian model 



We consider a chain of particles for simplicity in a one-dimensional lattice. In its most 
general form, the model is described by the following Hamiltonian 

N 2 

h(x, p ) = J2^ + a) (l) 

1=1 

where A = (Ai(x), A 2 (x), . . . ,A N (x)) is a generalized 'vector' potential in R N and x = 
(xi,x 2 , • • • , xjy), p = (pi,P2, ■ ■ ■ ,Pn) denote the particles position and momentum respec- 
tively. The Hamiltonian equations of motion read 

it = P - + A 

They are more transparent in Newtonian formulation 

dxi 

—r~ = 

dt 

dvi A 

dT = E^ + ir ( 3 ) 



where 



_ gA(x) _ &4 3 -(s) 

is an antisymmetric matrix containing in the entries the 'magnetic fields' acting on the 
particles. 



2.1 Conservation Laws 

Conservation of Energy: Even if the forces depend on velocities and positions, the model 
always conserves the total (kinetic) energy if the fields are time-independent due to the 
antisymmetry of the matrix B: 

i ( E = E B *w> = ( 5 ) 



Conservation of Momentum: Additional conserved quantities can be imposed by a suitable 
choice of the magnetic fields. For example, if we choose the Ai(x) such that they are 
left invariant by the simultaneous translations Xi — > Xt + 5, then the quantity J2iPi * s 
conserved. If in addition we require that ^2 i A i (x) = then this also implies the 
conservation of Yl,i x i- 



3 Time-independent Hamiltonian, continuous time 



The simplest way to realize this is the following: for a system with only energy conserva- 
tion, we may put Ai = fi(xi, Xi+i) for some functions For periodic boundary conditions, 
one identifies the sites modulo N. Then, the matrix B takes the nearest-neighbour form: 
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An open chain is obtained putting Bn = 0. A momentum-conserving model with next- 
to-nearest neighbor interaction can be obtained with: 



Ai — Cj+i — Cj 

Cj = fi(%i+l ~ x i) 



(7) 



for any functions /j. Again, one can take periodic boundary conditions (and then the 
indices are interpreted as modulo N), or consider an open chain, in which case one has 
to make C\ — 0. 

In both cases the end sites may be coupled to Langevin baths by modifying the equa- 
tions of motion of p± and pn adding a noise and a friction term, satisfying the fluctuation- 
dissipation relations for temperatures Ti and T N , respectively. 



4 Discrete time 



One can also consider discrete-time versions of these models, yielding a system of cou- 
pled maps. These have the advantage of being particularly easy to simulate, and, more 
importantly, they are analytically solvable in the high energy limit. As we shall see, the 
price we shall pay is that the map will not be simplectic [H], but it does preserve the 
phase-space volume and the Gibbs measure. 

We first describe a version without momentum conservation in detail, and then de- 
scribe more briefly the momentum-conserving case. We consider impulsive magnetic fields 
which act periodically and are of the form: 

Bij(x, t) = B% en (x)K even (t) + B°f{x)K oM {t) (8) 

where 

+00 

K even (t) = A(t-nr) 

n=— 00 
+00 

K odd (t) = Yl A(* - (n + l/2)r) (9) 

n=— 00 

where A(u) is a square pulse of intensity k and duration 1/k, and we consider k — > 00. 
During the short time when the fields are acting, the particles move a negligible amount, 
but their velocities rotate under their action. On the contrary, between kicks the velocities 
are constant and the motion is uniform. Without loss of generality we can set r = 1 
because changing the time interval between kicks amounts to rescaling the velocities. In 
the spirit of Refs. [23 EM EH] , we consider a chain with even N and choose the BfV en and 
B° dd such that the impulsive fields couple sites pairwise as follows: 

• B° dd is of the form © with B 2 = B 4 = B 6 = ... = 0, and B x = G{x x ,x 2 ), 
B 3 = G(x 3 ,x 4 ), etc. 

• B?? en is of the form © with B x = B 3 = B 5 = ... = 0, and B 2 = G(x 2 ,x 3 ), 
B4 = G(x4,X5), ... etc. 
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A periodic chain might be used, or, as we shall do later, one can connect the end sites 
to thermal baths during a half-cycle between kicks. The definition is completed by saying 
that whenever a site is connected to a heat bath during half a cycle, the interaction is so 
strong as to completely equilibrate the site at the bath's temperature. 

As we mentioned above, the map as we define it is not symplectic, even though it 
might seem as the composition of Hamiltonian steps. The reason for this is subtle: just 
before the end of the period during which the magnetic is on, the velocity is v i = p, + Ai, 
and just after it is turned off it is v[ = Pi. If we wish to impose that Vi = v[ (which we need 
for energy conservation), we see that this implies a jump in the momenta of magnitude 
Ai(x). This is a transformation that is not symplectic unless |^ = as is easy to 
check. An alternative way to see this issue is to consider the map as resulting from a 
time-dependent field, and then adding a (nonconservative!) force that exactly cancels ^ i 
in (jSJ), so that energy is conserved. Note, however, that the map conserves the volume for 
all trajectories, and leaves the Gibbs measure invariant, even if there is a potential V(q) 
acting during the intervals between kicks. 

Because in this time-dependent version we have given up symplecticity, we might just 
as well give up the requirement that the Bij = Bji derive from a 'vector potential' Ai. 
This will neither alter the equilibrium nor the volume conservation properties of the map, 
but gives us the freedom to use simpler functional forms for the B^. We tried different 
possibilities for the x-dependence of the magnetic field acting on each couple of sites, 
obtaining similar results. We present here those corresponding to the choice: 

G(xi, x i+1 ) = K[(xi + x i+ i) - 2tt] (10) 

With this choice the magnetic fields take values in the interval [— 2K-K, 2Ktv], where K is 
a constant parameter. 

By integrating the equations of motion (JHJ) between two successive kicks, and having 
imposed the condition of continuity of the velocities described above, we arrive at the 
following discrete time dynamics: denote by Ri(t) the 2x2 rotation matrix 

m - ( cW m) ) (ID 

9 



where 

c{Bi) = cos(Bi(t)) 
5(5,) = sm(Bi(t)) 

and consider the matrices C(t) and D(t) composed of a series of two by two blocks 
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Given x(t),v(t), the position-velocity vector at time t, its evolution will be 



x(t + l/2) = x(t) +v(t) 
v(t+l/2) = C(t)-v{t) 



(15) 



x(t + l)=x(t + l/2)+v(t + l/2) 
v(t + l)=D(t)-v(t+l/2)+£(t) 



(16) 



where £(t) is a iV- dimensional vector which is introduced to model interaction with heat 
baths 











(17) 



V mo j 

At each time the non-zero components £i(t) and £at(£) are i.i.d. random variables with 
Gaussian distribution 



p(v) 



i 



y 

-e 2^ 



(18) 



V2n<7 2 

where the variances of £i(t) and ^Ar(t) are respectively Ti and T N . 

All odd sites interact with the site to their right during half a cycle, and with the site 
to their left during the other half cycle. The sites at the ends interact with the baths 
during the corresponding half-cycles, and during that time they completely thermalize. 
In order to keep things simple we restrict the configuration space of each particle to the 
torus (x G [0, 27r]), so that all the Xi(t) are understood modulo 2n. 

As usual with these maps, if we consider the limit of small velocities and weak kicks, 
we recover a continuous evolution. 
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4.1 Chaoticity properties of the map. 

In order to see what chaoticity properties to expect from a chain, we begin in this section 
with the study of the dynamical properties of the elementary map with only two sites 
(x\, X2) = (x, y). This a particle moving in the x—y plane (actually we restrict ourselves to 
the torus T 2 ) under the action of an impulsive field along the z direction whose amplitude 
depends on the particle's position itself. Between kicks the motion is free; at each kick 
the velocity vector is rotated by an angle which is exactly given by the magnetic field, 
evaluated at the point where the kick takes place. Since the dynamics conserves the 
energy the accessible phase space is 3-dimensional. Denoting by v the modulus of the 
velocity (v = ^Jv\ + v 2 ) and by (3 its angle w.r.t. the x direction in the plane, we have 

x(t + l) = {x(t) + v cos((3(t))} mod (2vr) 
y(t + l) = {y(t) + vsm(f3(t))} mod (2tt) 

(3(t + l) = p(t)+K[(x(t + l)+y(t + l))-2ir] (19) 

The two-site map has two control parameters: the modulus of the velocity v (a constant 
of motion) and the amplitude of the field K. To illustrate the effect of varying them we 
have studied the Poincare section. In Fig. (|4.1|) we plot the surface section obtained by 
using the plane y = for K — 1 and several initial conditions. We see that as v tends to 
zero the trajectories are regular. As v increases, a smooth transition to a chaotic behavior 
is observed. If one considers a chain of more than two sites of coupled maps of this kind, 
one can expect the property of chaoticity to be stronger for larger system sizes. 

4.2 Numerical Analysis of the Fourier Law 

To compute the heat conductivity one has two possible procedures: the first is the direct 
non-equilibrium measure with two reservoirs at different temperatures, taking the ratio 
between the time averaged flux and the temperature gradient. If a linear response regime is 
applicable, one can also use the Green-Kubo formula, which enables to calculate transport 
coefficients as integrals of autocorrelation function in equilibrium states. The thermal 



12 




Figure 1: Poincare sections for the map (|19|) . Here K = 1 and v = 0.01 (top left), v = 0.6 (top 
right), v = 1. (bottom left), u = 5. (bottom right) 
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conductivity is then given by 



N 

Ji2 



oo 



k(T) 



< J(t)J(0) >GB (1 



2 



(20) 



where J(t) is the flux per particle 



(21) 



The local heat flux </$(£) is given by the change in energy for the couple (i, i + 1) before 
and after a kick. The angular brackets denote here an equilibrium average and the factor 
with the Kronecker delta function arise from the discreteness of time. 

We first performed the numerical simulation with two heat baths at temperatures 
Tl = 0.08 and Tr = 0.12. After a transient of 10 7 steps, we checked that a linear 
temperature profile is established by measuring the temporal average of twice the kinetic 
energy at each site (see Fig. (|4.2j) -left). We also kept track of the time averaged flux until 
it converged to its stationary value. We used at least 10 8 steps to check that fluctuations 
are less then a few percent. The measured heat conductivity is reported in Fig. (|4.2Jl -right. 
for different chain lengths from N = 8 to N = 2048. One can see that the conductivity 
has some fluctuation in the interval [0.37, 0.38] with an overall trend to stabilize, which 
suggests the approach to a finite value in the thermodynamic limit. We ran the dynamics 
for a large number of different initial random conditions (at least 10 3 ) and computed 
the heat conductivity by the Green-Kubo formula. We chose e = E/N = 0.05 which, 
according to the virial theorem, gives a kinetic temperature T = 0.1, the average value of 
the temperature of the two baths. The best fit of the data shows a convincing evidence 
that finite size effect are of the form 0(1/ N) with an asymptotic value = 0.376. 

4.3 High temperature limit of the discrete model 

In the high temperature limit the velocities Xi are large so that between two consecutive 
kicks new positions are translated by a large amount. Since the spatial coordinates are 
taken modulo 2tt, the sequence of positions constitutes a (quasi) random number gener- 
ator. Hence, we assume that in this limit, the position of the particle when the fields 
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Figure 2: Results of simulations for K = 1. Left: temperature profile with two heat baths for 
iV = 16 (circles) and ./V = 2048 (full line). Right: Heat conductivity versus the chain length 
N obtained throught non equilibrium simulations (squares) and GK formula (triangles). The 
continous line represent the best fit with a function koo + a/N. 

are flashed can be taken as uniformly randomly distributed. Because the magnetic fields 
are functions of the positions, this means that the fields themselves are random, with a 
distribution that is easily derived from the position dependence of the fields. All in all, 
what we have is that the velocity vector turns at each kick by a random amount, the 
probability distribution of which is known. Let us consider the change in energy for the 
couple + 1) before (e* = \vf) and after (e£ = \v' 2 ) a kick. 

e! { = c 2 {B i )e i + s 2 {B i )e i+1 + 2s{B i )c{B i )^^ 
e' i+1 = s 2 (B i )e i + c\B i )e i+l -2s(B i )c(B i )^^- l (22) 

where, of course, e[ + e' i+l = + e,i + \. Suppose now that Xi, a^+i are well approximated 
by independent uniform random process on the interval [0, 2tt]. Then the magnetic field 
(jlUj) will be a random variable and, with this particular choice for G, its probability 
distribution in the interval [—2Ktt, 2Ktc] will be 

f 5F- (! + x e \-2Kir, 0] 

p(B) = { 2KnK 2K * > L J (23) 



500 1000 1500 2000 

N 
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The mean energy of a couple of particles will be redistributed according to the rule: 

< e • > = x K < et > + (1 - x K ) < e i+ i > 
< e' i+1 > = (1 - x K ) < a > + x K < e i+1 > (24) 



since 



-2K* 2 1 V 2Kn 




< > = I coAB),,{B)dB = - ( I + ( : o ;.„ • ) ) -r,v 




V 2iT7r 

<s(S f )c(Si)> = (25) 

Here expectation values are taken over an ensemble of systems having the same velocities 
but random positions. If we choose K to be an integer or a half-integer, the dynamical 
rule (|24|) becomes [only as far as the means are concerned) the one of the stochastic model 
introduced by Kipnis, Marchioro and Presutti [23 E01 EI], where the total energy of a 
particle pair is equally redistributed between the two (x = 1/2). 

Let us give here a short computation of the heat conductivity in the case K = 1. The 
idea is to use first self-averaging with respect to the noise and then stationarity. The local 
'temperature' is defined as twice the kinetic energy Tj(i) = < vf{t) >. Using (|15 p -()16j ) and 
imposing stationarity Tj(i) = Tj(t + 1) = Tj one obtains the following recursion relation: 

Ti = T L 

-T<_ 2 + 27} - T i+2 = i = 2, 4, 6,..., iV -2 

T i+ i = Ti i = 0,2, 4, 6,..., N 

T N = T R (26) 

This can be easily checked starting from a configuration (Tl, T 2 ,T 2 ,T4,T4, T/v_ 2 , T/v_ 2 , Tr), 
making it evolve through the two kicks, and demanding that the same configuration is 
recovered. From the solution of the system (J26j) . a linear temperature profile is obtained 
(apart from the fact that the sites are in pairs): 

= T L + (Tj? ~ Tl) z (27) 
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It is easy to write a time-dependent version of (|26|). and check that the convergence to 
this profile is rapid, independently of the length of the chain, as it is essentially a diffusion 
equation. 

Between t and t + 1/2 we consider the local flux of the odd sites, while between t + 1/2 
and t + 1 we consider the local flux on the even sites. Self-averaging plus stationary gives 

Jo = \{T L + TJ-T L 

3i = -(Ti + T^-Ti z = l,2,...,iV-l 
Jn = \{T R + T N )-T N (28) 



The solution is an average site-independent local flux 

~N 



* = V- 1 ^ (29) 



and a spatial average flux 

N 

J 



N 

1 x 1 Tr — Tl 



iV^ J 2 iV 

i=0 

Putting together ()27j) and 1)30 j) we verify Fourier law J = kVT with an heat conductivity 
k = 1/2. 

In Fig. ()4.3|) we report the result of microcanonical simulations at different temperatures. 
One can see that for temperatures T > 10 the heat conductivity is basically constant 
and its numerical value coincides with the value we have just calculated in the high 
temperature regime. 

4.4 Momentum-conserving model 

We now run the arguments for a momentum-conserving model. Consider a coupling 
between three sites which is given by a magnetic field in the direction (1,1,1). The 
analogue of (|TT|) for this case is: 

' 1 + 20(73,) l-c(Bi) + V3s(Bi) 1 - c(A) - V5s(Bi) ^ 

l-c(Bi)-VZs(Bi) 1 + 20(73,) 1 - c(73,) + s(Bi) 

K 1-0(73,) + ^ 5(73,) 1-0(73;)- ^3 5(73;) 1 + 20(73,) ) 

(31) 
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Figure 3: Thermal conductivity k versus the temperature (note the log-log scale). The dashed 
line corresponds to the value k = 1/2. Results are plotted for three different N values (see leg- 
end). Here only a single trajectory has been used in the application of the Green-Kubo formula. 
However the good overlap between the three curves indicates reliability of the computation. 
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where now 

c{Bt) = cos(V3Bi(t)) 
s(Bi) = sm(y/3Bi(t)) 

A simple choice for B^ = G(xj, 0^+2) is: 

G(xi, x i+ i, x i+2 ) = K[(xi + x i+ i + x i+2 ) - 3ir] 



(32) 



(33) 



Clearly, the velocity in the (1,1,1) direction (the sum of the three velocities) is not affected. 
If we now alternate kicks which couple sites (1,2,3), (4,5,6), (7,8,9) etc, with kicks 
which couple sites (2,3,4), (5,6,7), (8,9,10), we obtain a map that at each kick 
conserves the total sum of the velocities (and, of course, the energy). 
We can repeat the argument in the previous section to obtain the high energy limit. We 
set the multiplicative constant K to 1/V3, for simplicity. If we consider that at each step 
the positions x i: x i+ i, x i+2 are uniform and independent random with the choice ()33j) we 
have that: 

(c(B t )) = (s(B t )) = (ciBiWBi)) = (34) 



(c 2 (B t )) = (f{B t )) = - 
and we see that the average energies evolve during a kick as: 



' (ei) ^ 
«+i> 

V ( e 'i+2) J 



/ 



\ 



1 1 
111 
1 1 1 



1 \ / (e.) N 
(e*+i) 



(35) 



(36) 



Again we can easily find an equation for the T», and check that the profile is linear, with 
the only difference that now the sites come in triplets: 



/ rri rri rri rri rri rri rri rri rri rri rri rri 

[J-L, J-L, J-3, J-3, J-3, J-6, J-6, J-6, ■ ■ ■ , N—3: J-N-3, J-N-3, J- J 



R 



(37) 



and the temperatures satisfy 



-Ti_ 3 + 2T t - T i+3 = % = 3, 6, 9, . . . , AT - 3 



(3? 
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Figure 4: Same as in Fig. (|4.H|) but for the momentum conserving case. The analytical high- 
temperature value is k = 1. 

From this we easily obtain a value for the thermal conductivity k — 1. 

The numerical simulations do agree with this value in the high temperature regime 
and give a finite heat conductivity at each temperature, see Fig. 14.41 

5 A stochastic minimal model 

As was done in the previous section for the discrete time model we study here a high 
temperature limit which becomes a continuous time model, by replacing the magnetic 
fields with suitable random processes. 

5.1 Fokker-Planck equation 

In the limit in which the field is very weak, and the velocities (and hence the energy) are 
very large, we get a velocity field that is randomly exchanging momentum, but in small 
amounts each step. The deterministic equations of motion © can then be replaced by a 
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system of Langevin equations with multiplicative noise 

^ = X>^ (39) 
with the Bij of the nearest-neighbour form 

By = B i (t)(5 iJ+ i-5 i+ i J -) (40) 
and Bj a white, Gaussian random variables with 1 

< Bj(t) > = < Bj{t)B k {t') >= 2 5(t - t') (41) 



Standard techniques [o9] yield the following Fokker-Planck equation for the evolution on 
the velocities probability distribution 



dP 



( 42 ) 



dt 

i=i 

where the operators L i)i+ i are 'angular momentum' operators corresponding to the Lapla- 
cian on the sphere: 

d d 

L Li+1 = Vi v i+1 — (43) 

ov l+1 dvi 

If we add to the bulk system the interaction with two heat baths connected to the first 
and last particle (respectively at temperatures T L and T R ) and represented as Ornstein- 
Uhlenbeck processes, we arrive at 

Qp N 

— = L h+i p + diKB) + T L d\P + d N (v N P) + T R d 2 N P (44) 
i=i 



This diffusion process has also been considered by Olla [40] , in combination with a deter- 
ministic harmonic part. 

One can also easily write a momentum-conserving variant. Indeed, considering the 
low-field version of the conserving map, we get: 

OP N 

= ~ ^^{-^M+l + -^i+l,i+2 + L i+2t i} 2 P (45) 

8=1 

where in a closed chain the indices are understood modulo N. 



lr This being a multiplicative process, the convention (Ito, Stratonovitch...) should be specified. In 
fact, the precise definition will be implicit in the Fokker-Planck equation (|42|l 143(1 
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5.2 Expectation values 

We would like to compute the energy and flux expectation values and also energy cor- 
relations. One may derive an equation of motion for these quantities directly from the 
Fokker-Planck equation. Multiplying for instance (pf4*|) by v 2 and integrating the resulting 
equations, i.e. 

| J v 2 P(v)d N v = J v 2 

+ j v 2 (d^Piv)) + T L d 2 P{v)) d N v 

, 2 (d N (v N P(v)) + T R d 2 N P(v)) d N v (46) 



by using repeatedly integration by parts (with Vidi = diVi — 1) we obtain 

^ < v\ > = -4 < v\ > +2 < v\ > +2T L 



— <v 2 i> = 2 < v\_ x > -4 < v\ > +2 < v 2 i+1 > for i = 2, . . . , N — 1 



d_ 

It 

J t <v N> = - 4 < V N > + 2 < V N-l > + 2T R (47) 

In the stationary state j| < vf >= so that the temperature profile is obtained as the 
solution on a linear system of equations. The solution is: 

T l = T L + (T ^~^ L) i (48) 

In the stationary state there is a net heat current J flowing through the lattice. This can 
be calculated by directly measuring the energy exchange with the two baths. The energy 
flux from the left reservoir to the first particle is 

<J 1 > = <v 2 >-T L (49) 

while the energy flux from the last particle to the right reservoir is 

<J N > = <v 2 N >-T R (50) 

In both cases, by using we obtain 

<J 1 > = <J N > = J= ^'J^ (51) 
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Putting together ()48|) and (joTj) we find that Fourier law holds for the stochastic model 
with a heat conductivity k = 1. 



5.3 The stationary measure 



We show here that the stationary measure in the presence of heat flow cannot be of either 
Gaussian or product form. First of all, the evolution operator is invariant with respect 
to change in sign of any single velocity, and hence P{yi . . .v^) = P(T%, ...,T N ), i.e. the 
velocities only appear as even powers. 
Proposing a Gaussian measure of the form: 



P(vi...v N ) ~ e ~^' A ^ 



(52) 



substituting into (jSJ), we get: 



n 2 



t jViVj 



ijViVj 



IE 



vf) VveR N 



Going to the diagonal basis A = Ai5ij, this equation implies (Ai — A\ 
proportional to the identity. This is impossible if there is heat flow. 
Let us see that a product measure: 



(53) 

0, i.e., A is 



TV 

P(vi ...V N ) = \\pi{Vi 

i=l 



(54) 



is in general not possible if there is heat flow. In the stationary state we have 

LppP 
P 



L M+l^(^)^+l(^+l) 







(55) 



Pi(Vi)p i+1 (Vi + i) 

Because of their different arguments, each term must vanish separately, thus: 

(£r?_ M + Ll i+1 )p i - 1 (v i -i)p i (vi)p i+1 (v i+1 ) = (56) 

This is an equation of the form [L\ + Ly)ip = (with L x , L y the SU(2) operators), which 
can only be satisfied by the quantum numbers (l,m) = (00), the spherically symmetric 
function. In our case this means that: 



Pi(x)p 2 (y)p3(z) = g(x 2 + y 2 + z 2 



(57) 
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which is a product only if we have an isotropic Gaussian, again impossible in the presence 
of heat flow. 

The equation of motion for n-point correlation functions of these models are linear and 
close within themselves. This is a strong suggestion that the models may be integrable, 
but we shall not persue this line here. 

6 Conclusions 

We have studied a family of models of heat transport whose high energy limit is ana- 
lytically solvable. The intermediate energies are easy to simulate, the approach to the 
asymptotic value can be easily tested. The momentum-conserving version of the model 
has finite conductivity, providing further confirmation of the fact that in itself momentum 
conservation does not imply anomalous conductivity. 

It would be nice in the future to investigate in the simple framework of "magnetic 
kicks" we have introduced a system which has also a potential energy. For example one 
could start from a system of coupled linear oscillators, establishing a connection with the 
work of Olla 1401. 



Another interesting problem is to solve the Fokker-Planck equation to find the sta- 
tionary distribution of our system. This could shead some light also on the stochastic 
model introduced by Presutti et al. [3U]. There again the stationary measure is a product 
measure only locally (i.e. in the approximation to first order in L^ 1 ). 
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